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P-CONSISTENCY  ^ 

BY 

Patrick  J.  Rabier  and  Werner  C.  Rheinboldt^ 

Abstract.  Part  I  of  this  paper  presented  a  theory  of  distribution  solutions  of  semilinear 
differential-algebraic  equations  (DAE’s).  In  particular,  it  was  shown  that  uniqueness  of  solu¬ 
tions  of  initial  value  problems  breaks  down  completely  in  the  class  of  discontinuous  solutions. 
Here  a  mathematical  procedure  is  introduced  for  selecting  physically  acceptable  solutions 
which  satisfy  some  new  consistency  condition  relative  to  admissible  perturbations  of  the  orig¬ 
inal  DAE.  Several  nonlinear  circuit  examples  are  given  to  support  the  theory. 


1.  Introduction. 

As  shown  in  I,  Section  3^  of  this  paper,  the  difficulties  associated  with  the  existence  of 
impasse  points  for  DAE’s  of  the  form 

(1.1)  A{i)x  =  G{t,x) 

can  be  partially  resolved  by  considering  distribution  solutions  instead  of  classical  ones. 
Unfortunately,  this  creates  the  new  problem  of  invalidating  the  uniqueness  of  solutions  to 
initial  value  problems  for  (1.1)  which,  in  general,  will  now  have  not  one  or  even  a  few 
solutions,  but  uncountably  many  of  them  (see  I,  Example  3.1).  Of  course,  this  mathe¬ 
matical  fact  is  hai-dly  acceptable  on  physical  grounds,  and  can  only  mean  that  most  of 
the  solutions,  optimally  all  but  one,  are  deprived  of  physical  significance  and  should  be 
discarded.  This  raises  the  question  about  suitable  criteria  for  such  a  selection. 

^The  work  was  supported  in  part  by  ONR-grant  N-00014-90-J-1025,  and  NSF-grant  CCR-9203488 
^Department  of  Mathematics  and  Statistics,  University  of  Pittsburgh,  Pittsburgh,  PA  15260 
^Throughout  the  presentation  this  notation  refers  to  the  corresponding  item  in  Part  I  of  this  article. 
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At  first  sight,  this  task  seems  to  require  the  specification  of  some  physical  setting  for 
(1.1)  and,  in  particular,  of  appropriate  characteristics  of  the  physically  acceptable  solu¬ 
tions.  Here,  we  shall  develop  the  thesis  that  a  purely  mathematical  selection  procedure 
can  be  introduced,  and  that  physical  data  are  needed  only  when  it  comes  to  the  actual 
interpretation  of  the  results.  We  begin  with  an  outline  of  this  idea. 

Sometimes  the  opinion  has  been  fostered  that  every  DAE  should  or  could  be  considered 
as  the  limiting  case  “c  =  0”  of  a  singularly  perturbed  ODE.  However,  physically  motivated 
perturbations  generally  lead  to  such  ODEs  only  after  considerable  manipulation  of  the 
equations  and  under  global  invertibility  assumptions  that  are  not  necessarily  satisfied  by,  or 
very  difficult  to  check  in,  nontrivial  concrete  problems.  In  particular,  the  simple  procedure 
of  changing  the  algebraic  constraints  in  (1.1)  into  differential  equations  weighed  by  the 
factor  e  usually  results  in  problems  totally  deprived  of  physical  reality,  and  hence  of  no 
help  for  a  better  understanding  of  the  case  e  =  0.  In  fact,  what  is  often  apparent  is  that 
the  DAE  (1.1)  does  represent  a  simplified  model  of  some  physical  system  in  which  small 
parameters  have  been  neglected,  but  that  the  perturbed  problems  remain  DAE’s  (and 
hence  may  be,  though  rarely  are,  explicit  ODE’s).  Another  difference  is  that  the  neglected 
parameters  form  a  vector  which  reduces  to  a  single  scalar  only  in  the  simplest  models. 

Accordingly,  we  shall  assume  that  (1.1)  arises  as  a  simplification  of  a  perturbed  DAE 
obtained  by  disregarding  some  physical  parameters  of  negligible  magnitude  represented  by 
a  vector  /i.  Typically,  the  exact  value  of  p  is  unknown,  but  the  //-dependent  perturbed 
problems,  henceforth  denoted  by  where  “P”  stands  for  “perturbation  ,  can  be  written 
down  explicitly  since  the  nature  of  the  neglected  quantities  is  clear  from  the  context.  In 
this  respect,  see  Remark  1.1  below. 

With  every  fixed  value  of  the  small  perturbation  parameter  //  we  shall  associate  an  open 
subset  W^'’^  of  the  set  of  consistent  points  of  (1.1)  (see  Part  I)  depending  only  upon 
the  direction  of  //  and  consisting  of  points  with  the  following  property:  If  (fo,a:o)  € 
then  every  solution  of  (1.1)  passing  through  xq  at  time  to  can  be  approximated  imiformly 
in  some  open  interval  about  to  by  solutions  of  the  perturbed  problem  as  €  — >  0’*’.  Thus, 

if  n  is  small  in  norm  in  the  first  place,  solutions  of  (1.1)  passing  through  points  of 
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are  already  closely  approximated  by  solutions  of  when  €  =  1;  that  is,  by  solutions  of 
The  set  will  be  characterized  in  terms  of  an  eigenvalue  sign  condition  which  also 
guarantees  the  aforementioned  (local)  approximation  property  via  the  Tikhonov-Levinson 
theorem  of  singular  perturbation  theory  and  its  generalization  to  DAE’s  in  the  dissertation 
[Ya].  Points  of  will  be  called  V ft- consistent  points. 

Our  criterion  for  the  selection  of  the  physically  relevant  solutions  of  (1.1)  is  based  upon 

(i)  the  reasonable  assumption  that  such  solutions  should  be  close  (at  least  as  distributions) 
to  solutions  of  Vfi  if  indeed  Vft  represents  a  realistic  small  perturbation  of  (1.1),  and  upon 

(ii)  the  postulate  that  physical  solutions  of  (1.1)  should  not  exhibit  discontinuities  (that  is, 
they  should  be  classical  ones)  unless  their  continuity  creates  a  dichotomy  with  (i)  above. 
This  premise  may  perhaps  be  disputed  in  some  problems,  but  we  are  not  aware  that  its 
plausibility  has  thus  far  been  challenged.  In  fact,  it  is  implicitly  made  in  all  the  traditional 
singular  perturbation  approaches. 

The  arguments  developed  so  far  provide  an  immediate  justification  to  our  procedure  for 
selecting  the  solutions  of  (1.1)  compatible  with  the  perturbed  problem  "P^.  By  definition, 
a  solution  of  (1.1)  will  be  called  -  consistent  if  {t,x{t))  €  for  all  t  in  the  time 

interval  of  interest,  and  if  no  jump  occurs  at  auiy  point  of 


It  is  important  to  notice  that  while  the  condition  that  no  jumps  occur  at  points  of 
eliminates  many  otherwise  possible  distribution  solutions,  it  does  not  rule  out  the 
possibility  of  jumps  altogether,  which  may  still  exist  when  a  solution  reaches  points  of 
the  bormdary  \  W^’^.  Furthermore,  a  jump  will  exist  if,  in  addition,  there  exists  a 

P^-consistent  point  where  it  terminates  (otherwise,  the  solution  is  not  extendable  to  larger 
time).  In  particular,  jumps  may  always  take  place  from  singularities  such  as  impasse  points, 
which  are  not  even  consistent  points  (i.e.  points  of  As  we  shall  see,  jumps  are  also 

possible  from  points  which  are  not  singular  points  of  (1.1)  and,  in  fact,  which  cannot  be 
characterized  independently  of  the  perturbation. 

Remark  1.1:  For  the  claim  that  the  solutions  of  (1.1)  selected  by  the  P^- consistency 
criterion  are  physically  meaningful,  it  is  of  some  importance  that  the  components  of  the 
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vector  account  for  all  the  physical  parameters  neglected  in  the  original  problem.  The 
difficulty  in  re-introducing  these  parameters  varies  with  the  nature  of  the  physical  prob¬ 
lem.  For  instance,  in  elementary  particle  mechanics,  masses  of  rods  (and  sometimes  of 
some  particles)  are  commonly  omitted  from  the  models.  In  LRC  networks,  phenomena  in¬ 
duced  by  the  wiring  are  usually  disregarded,  and  they  can  be  accounted  for  by  introducing 
parasitic  series-inductainces  and  shunt-capacitors.  □ 

In  the  above  approach  the  only  obstacle  to  really  singling  out  the  perturbation  of  (1.1) 
with  actual  significance  is  the  fact  that  the  exact  value  of  the  parameter  //  is  not  known  in 
general.  On  the  other  hand,  our  examples  indicate  clearly  that,  qualitatively,  this  value  is 
relatively  unimportant,  and  that,  in  general,  the  solutions  are  observed  to  show  only  very 
few  distinct  behavior  characteristics  as  is  varied.  This  suggests  that  it  may  be  possible  to 
classify  perturbed  DAE  problems  in  terms  of  inequivalent  normal  forms,  in  the  very  same 
spirit  as  has  been  done  (to  some  degree)  for  perturbed  bifurcation  problems.  A  significant 
added  difficulty  is  the  necessarily  global  aspect  of  such  perturbed  DAE  problems  that  jump 
phenomena  make  no  sense  in  a  framework  of  local  analysis.  These  questions  will  not  be 
considered  further  here. 

Naturally,  the  idea  of  relating  the  solutions  of  DAE’s  to  those  of  perturbed  problems 
is  not  new.  The  eigenvalue  condition  central  to  this  Part  11  has  the  visible  appearance 
of  a  stability  condition  for  equilibria.  In  the  work  on  nonlinear  circuits  by  Green  and 
Willson  [GW]  stability  is  defined  only  for  stationary  points  (operating  points);  that  is, 
for  solutions  of  the  dc-circuit  obtained  by  disregarding  all  dynamic  elements.  Then,  m 
essence,  an  operating  point  is  stable  if  under  aU  possible  perturbations  by  means  of  parasitic 
capacitors  and  inductances  the  corresponding  stationary  point  of  the  resulting  dynamic 
circuit  is  also  stable. 

On  the  other  hand,  Sastry  and  Desoer,  [SDe],  consider  perturbations  of  a  given  circuit 
which  produce  a  singularly  perturbed  ODE  involving  a  scalar  perturbation  parameter. 
Then  stabihty  of  a  solution  of  the  unperturbed  problem  can  be  defined  by  means  of  the 
conditions  of  the  standard  Tikhonov- Levinson  theorem,  which  turn  out  to  be  closely  related 
to  those  in  this  paper.  One  major  difference  here  is  that  we  allow  perturbations  to  DAE 
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form  and  involving  a  vector  of  small  parameters.  Moreover,  in  our  theory  we  consider, 
from  the  start,  distribution  solutions  while,  in  [SDe],  discontinuous  solutions  of  (1.1)  have 
to  be  defined  as  the  limit  of  continuous  solutions  of  the  perturbed  problems.  Thus,  in 
our  approach,  the  solutions  of  (1.1)  are  defined  independently  of  the  perturbations  which 
become  involved  only  later,  in  the  selection  process.  A  major  technical  advantage  is  that 
we  do  not  require  any  assumption  about  the  behavior  of  the  solutions  of  the  perturbed 
problems,  and  thus  can  bypass  serious  difficulties  in  the  vicinity  of  jumps.  As  a  result, 
we  are  able  to  handle  situations  both  frequently  encountered  and  not  covered  by  the 
analysis  of  [SDe].  On  the  other  hand,  the  issue  whether  our  'P^-consistent  solutions  can 
be  approximated  in  the  sense  of  distributions  (and  hence  globally)  by  the  solutions  of  the 
problems  as  e  — >•  0+  is  an  important  question  as  yet  imresolved  by  the  segment  of  the 
theory  presented  in  this  paper. 

Section  2  provides  the  technical  details  of  the  preceding  discussion.  The  relationship  of 
P-consistency;  that  is,  P^  -  consistency  for  some  with  the  Tikhonov-Levinson  theorem 
is  examined  in  Section  3.  Our  aim  there  is  to  justify  the  choice  of  our  conditions  by 
showing  how  they  relate  to  this  classical  result  of  singular  perturbation  theory.  It  is  a 
significant  advantage  of  our  method  that  the  extremely  cumbersome  technicalities  involved 
in  this  process  never  have  to  be  worked  out  in  any  concrete  example.  In  contrast,  any 
approach  based  upon  eventual  reduction  to  a  standard  perturbation  problem  for  ODEs, 
when  possible  at  all,  would  have  to  display  similar  technicalities  in  every  case. 

In  Section  4,  several  examples  from  electrical  circuit  theory  are  given.  As  already  men¬ 
tioned,  in  such  problems  the  natural  perturbations  are  represented  by  parasitic  inductances 
and  capacitors,  which  turns  out  to  be  in  concordance  with  Smale’s  proposal  in  [Smj. 

The  notation  and  terminology  are  those  of  Part  I,  to  which  the  reader  is  referred  re¬ 
garding  all  symbols  and  concepts  not  defined  here. 

2.  P-consistency. 

As  noted  in  the  Introduction,  our  paradigm  is  that,  in  practice,  many  DAE’s  (1.1)  arise 
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“naturally”  as  limiting  cases  of  perturbed  DAE’s  of  the  form 

(2.1)  A{ti,t)x  = 

such  that  the  DAE  (1.1)  is  recovered  for  /i  =  0.  If  so,  A{t)  and  G{t,x)  in  (1.1)  become 
A(0,f)  and  G(0,f,x),  respectively.  Above,  by  “naturally”  we  mean  that  the  vector  6  R 
represents  a  collection  (mi  ,  -  •  •  ,  /xjt)  of  negligible  physical  constants.  An  apparently  innocu¬ 
ous  fact  will  turn  out  to  have  a  crucial  importance  for  the  justification  of  our  arguments: 
In  most,  if  not  all  cases,  physical  constants  have  meaning  only  for  positive  values  (after 
possible  rescaling  and  sign  changes).  As  a  result,  the  perturbations  (2.1)  involve  only 
vectors  fi  with  positive  components.  Accordingly,  we  shall  set 

(2.2)  =  {/i  €  R*  :  IImII  <}? 

where  =  (0,  oo)  and  k  is  a  nonnegative  integer  (when  k  =  0,  no  physical  constant 
has  been  neglected  in  (1.1);  this  case  does  fall  within  the  framework  of  this  paper  and 
is  discussed  in  Remark  2.4  later).  In  (2.2),  the  normalization  \\n\\  <  1  could  be  replaced 
by  ||/i||  <  V  where  77  >  0  is  arbitrary.  The  choice  r?  =  1  is  made  here  for  convenience, 
and  mainly  to  avoid  introducing  a  parameter  rj  which  would  never  be  involved  in  further 
considerations. 

Remark  2.1:  In  (2.2),  the  choice  is  justified  because  it  is  assumed  that  /x 

represents  only  physical  parameters  which  are  ignored  in  (1.1)  (and  hence  set  to  0  there). 
This  is  an  important  point,  for  if  /x  includes  components  representing  perturbations  of 
physical  constants  already  present  in  (1.1),  then  negative  as  well  as  positive  values  of  these 
components  are  perfectly  admissible.  There  is  no  difficulty  to  incorporate  perturbations  of 
already  existing  physical  parameters,  thus  letting  (2.1)  depend  upon  an  additional  vector  u. 
Under  conditions  frequently  met  in  applications,  doing  so  does  not  modify  the  final  results, 
and  the  behavior  of  interest  continues  to  be  solely  the  consequence  of  those  parameters 
that  are  not  already  embodied  in  (1.1).  Thus,  to  keep  technical  matters  as  simple  as 
possible  and  not  to  obscure  the  issue,  this  otherwise  legitimate  generalization  will  not  be 
considered  here.  □ 
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Standard  perturbation  theory  uses  fc  =  1  in  (2.2)  and,  hence,  a  positive,  real  parameter 
fi  (usually  called  e).  Moreover,  it  requires  that  (2.1)  is  an  ODE;  that  is,  A{fJ.,t)  is  invertible 
for  any  vector  /i  with  positive,  and  sufficiently  small  components. 

For  technical  reasons,  we  assume  that  the  mappings 

(2.3)  eC+x  €  £(K'‘) 

and 

(2.4)  G  X  J  X  R"  I — >  G  R”, 

are  of  class  C°°  (sufficient  smoothness  would  suffice)  and  that 

(2.5)  rank  A(/i,t)  =  p  >  r  =  rank  A(0,t),  V(p,t)  G  C+  x  J 
and 


(2.6)  ker  A(p,i)  C  kerA(0,t),  V(p,0  G  C+  x  J. 


Thus  p  is  independent  of  (p,f),  while  the  independence  of  r  from  i  was  already  assumed 
in  Part  I.  The  condition  (2.6)  can  be  weakened;  see  Remark  2.6  below. 

Evidently,  (2.1)  is  an  ODE  for  p  G  C+  if  and  only  if  p  =  n  in  (2.5),  so  that  (2.6)  becomes 
vacuous.  For  future  use,  we  note  that  (2.6)  implies  that 


A 

ker  — (0,t) 


V(p,t)  X  J, 


Indeed,  if  f  €  J  is  fixed  and  (ej)  is  a  sequence  of  positive  numbers  tending  to  0,  and 
hj^i...,hj^n-p  denotes  an  orthonormal  basis  of  ker  A(ejp,<),  then  we  have  hj^i  G  ker  A(0,  t) 


and  hence 


(2.8)  ~{A{ejp,i)  -  A{0,t))hj^i  =  0,  l<i<n-p,j>0. 

After  extracting  a  subsequence,  we  may  assume  that  lim  hj^i  =  hj,  1  <  i  <  n  —  p,  so  that 
the  vectors  hi,--  -  ,/in-p  form  an  orthonormal  subset  of  ker  A(0,t).  Now,  by  taking  the 
limit  in  (2.8),  we  find  that  {{dA/dp){0,t)p)hi  =  0,  l<i<n  —  p,  and  (2.7)  follows. 
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As  in  Part  I,  we  focus  only  on  the  case  when,  for  ^  =  0,  (2.1)  (that  is,  (1-1))  is  a 
geometrically  nonsingular  DAE  of  index  one.  Let  (t,  a:)  G  be  a  consistent  point  for 
(1.1),  so  that  €  GL{kevA{0,t),Z)  by  I,  Proposition  2.1,  where 

Q{t)  is  the  projection  used  in  I,  Section  2  and  corresponding  to  A{t)  =  A(0,t).  (In  other 
words,  Q(t)  projects  onto  the  fixed  complement  Z  of  rge  A(0,t),t  G  J).  For  notational 
convenience,  we  set 


(2.9)  A(/i,t,x)  = 


Q(i)D,C?(0,t,x)i, 


Q(f)(~(0,%)€£(kerA(0,t)). 


Definition  2.1.  For  fixed  n  gC+  the  point  {t,x)  G  is  consistent  with  the  perturbation 
(2.1),  or  Vfi- consistent  for  short,  if  the  operator  A{p,t,x)  in  (2.9)  has  p  —  r  eigenvalues 
(counting  multiplicities)  with  strictly  negative  real  part.  The  set  of  Vf^-consistent  points 
will  be  denoted  by  . 

Occasionally,  we  shall  refer  to  a  point  (t,x)  being  P-consistent  if  it  is  P^-consistent 
for  some  p  G  "P.  This  terminology  is  compatible  with  calling  (2.1)  for  fixed  p  €  C+  the 
“perturbed  problem  and  calling  P  the  collection  of  perturbed  problems  P^  for  p  G  C+. 
The  notation  also  reflects  the  obvious  fact  that,  not  only  what  value  the  parameter  p  has, 
but  also,  and  primarily,  how  p  enters  the  system  (2.1)  is  important  here. 

By  (2.7),  the  operator  A(p,t,x)  in  (2.9)  has  0  as  an  eigenvalue  with  geometric  multi¬ 
plicity  at  least  equal  to  n  —  p.  Hence,  Definition  2.1  applies  if  and  only  if  0  is  an  eigenvalue 
of  algebraic  and  geometric  multiplicity  exactly  equal  to  n  -  p  while  its  other  eigenvalues 
have  negative  real  parts.  Note  also  that  the  condition  of  Definition  2.1  is  valid  if  and  only 
if  it  holds  with  p  replaced  by  any  positive  scalar  multiple.  On  the  other  hand,  because 
of  A(-p,t,x)  =  -A(p,t,x),  Definition  2.1  would  be  meaningless  if  p  and  -p  could  both 
represent  admissible  perturbations,  unless  p  =  r  and  hence  A  =  0  in  (2.9).  This  is  ruled 
out  by  the  condition  p  €  K  which,  in  turn,  was  earlier  dictated  by  the  positivity  of  physical 
constants. 

Remark  2.2;  For  given  u  G  kerA(0,f)  we  see  that  v  A(p,t,x)  u  G  kerA(0,t)  is 


characterized  by  the  condition 

€  rge  A(t), 

which  is  clearly  equivalent  to  (2.9).  This  shows  that  and  hence  its  eigenvalues 

and  their  multiplicities,  are  independent  of  the  choice  of  Q{t)-  □ 

Proposition  2.1.  For  p  €  K,  the  set  W^’>^  is  open  in  W^. 

Proof.  Suppose  this  does  not  hold.  Then,  there  exists  a  r,)  e  and  a  sequence 

itj,xj)  e  W  with  lim(fj,rj)  =  (<*,2:,).  Since  is  not  F^-consistent, 

the  operator  A(fjL,tj,Xj)  has  at  least  n-  p  +  l  eigenvalues  (counting  multiplicities)  with 
nonnegative  real  part.  In  essence,  a  contradiction  is  reached  by  using  the  continuity  of  the 
eigenvalues  with  respect  to  the  operator  (see  [K]),  which  implies  that  A(/i,<*,x*)  has  at 
least  n  -  p  +  l  eigenvalues  with  nonnegative  real  parts,  whence  (t,,  x*)  ^ 

The  above  argument  requires  the  domain  of  the  operators  to  be  the  same,  which  here 
is  not  the  case  since  kerA(0,t)  varies  with  t.  This  apparent  difficulty  is  overcome  by 
considering  a  one-to-one  mapping  U{t)  €  £(R”“'^,1R")  with  range  ker  A(0,  f)  which  de¬ 
pends  continuously  upon  t  near  f,.  Existence  of  U  is  easily  derived  from  the  constancy  of 
dimker  A(0,f)(=  n  -r).  The  eigenvalues  of  A(/i,f,  x)  for  (f,x)  near  (f*,x,)  are  also  those 
of 

[Q{t)D,G{0Xx)U{t)]-^Q{t)(^^(0,t,^)y{t)  e  £(R"-0, 
and  the  continuity  argument  may  now  be  used  with  these  operators.  This  completes  the 
proof.  □ 

The  notion  of  -consistent  point  leads  at  once  to  that  of  -consistent  solution: 

Definition  2.3.  Let  /i  6  £+  he  given.  The  function  x  ;  J  — >  R”  is  a  -  consistent 
solution  of  (1.1)  if 

(i)  (f,x(t))  €  IV"  (^x(t))  e  iv"'". 

(ii)  For  to  ^  J  and  {to,x{io))  G  IV"''',  there  is  a  T  >  0  such  that  coincides  with 

the  unique  classical  solution  of  (1.1)  passing  through  x(to)  for  t  =  to  (see  I,  Theorem  2.1). 
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The  union  of  all  the  V ^-consistent  solutions,  as  fi  €  is  varied,  is  called  the  set  of 
V-consistent  solutions. 

By  Definition  2.2,  a  P^-admissible  solution  x  of  (1.1)  cannot  jump  from  points  of 
Nevertheless,  it  may  have  a  discontinuity  at  a  point  (to, 2^)  €  \  if  e 

for  t  <  to  and  lim  a:(t)  =  a:o  •  This  may  happen  when  x^  €  (see  Remark  2.3  later 

inl¬ 
and  Example  3  of  Section  4),  or  when  Xq  is  an  accessible  impasse  point  (see  I,  Theorems  3.1 

and  3.2)  provided  that  there  is  a  (to,xJ)  €  such  that  —  Xq  €  kerj4(0,fo).  Jumps 
may  also  occur  at  other  singularities,  but  then  it  must  be  proved  that  the  discontinuous 
function  so  obtained  solves  (1.1)  in  the  sense  of  distributions.  This  amounts  to  justifying 
the  validity  of  formula  I  (3.1)  which,  as  noted  in  Part  I,  should  not  be  taken  for  granted. 

Because  (t,  x{t))  must  lie  entirely  in  when  x  is  a  "P^-consistent  solution  of  (1.1),  the 
occurrence  of  discontinuities  is  considerably  restricted.  However,  this  does  not  necessarily 
reinstate  the  uniqueness  of  solutions  of  the  initial  value  problems.  In  fact,  for  a  point 
(io,Xg)  e  \  the  “new”  starting  point  need  not  be  uniquely  determined 

in  As  noted  at  the  end  of  Part  I,  there  will  be,  in  general,  only  a  finite,  or  at 

most  countable,  number  of  possible  choices  for  x^ .  Of  course,  it  may  also  happen  that  no 
choice  for  x+  is  available,  and  hence  the  solution  x(t)  stops  at  x“  and  cannot  be  continued 
further. 

We  note  also  that  because  of  the  openness  of  the  set  (Proposition  2.1)  along  with 

condition  (ii)  of  Definition  2.2,  -  consistent  solutions  of  (1.1)  can  exhibit  jumps  only  at 

isolated  (and  hence  discrete)  time  values.  Thus,  altogether,  we  find  that,  in  general,  there 
are  at  most  countably  many  (and  often  only  finitely  many)  -consistent  solutions  of  (1.1) 
with  prescribed  initial  value  defined  in  the  open  interval  J . 

Remark  2.3:  From  the  above  comments,  P^  -  consistent  solutions  of  (1.1)  are  solutions 
which  are  classical  ones  unless  they  have  no  other  choice.  Reasons  for  exhibiting  dis¬ 
continuities  are  either  intrinsic  (singularities)  or  perturbation-dependent  (two  conjugate 
eigenvalues  of  the  operator  A(^,f,x)  in  (2.9)  cross  the  imaginary  axis  in  the  “wrong”  di¬ 
rection;  that  is,  pass  from  negative  to  positive  real  part).  In  both  cases,  jumps  occur  from 
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the  impossibility  of  approximating  would-be  classical  solutions  of  (1.1)  by  classical  solu¬ 
tions  of  Viti  (i-e.  (2.1)  with  /i  replaced  by  eji)  as  e  ^  0"^  in  the  viciiuty  of  the  point  (to,  ) 
where  the  jump  originates.  The  relationship  between  T’^-consistency  and  approximation 
will  be  investigated  in  detail  in  the  next  section.  For  the  time  being,  note  only  that  a 
solution  x(t)  of  (1.1)  can  certainly  not  be  approximated  by  any  continuous  fvmction  in  any 
open  interval  centered  at  to  when  Urn  x{t)  =  and  (to,a:o  )  is  an  (accessible)  impasse 
point  of  (1.1).  Indeed,  in  this  case,  there  is  no  continuous  extension  of  x{t)  as  a  classical 
solution  of  (1.1)  beyond  t  =  to.  □ 

Remark  2.4:  It  is  not  without  interest  to  apply  this  theory  to  the  case  when  p  =  r  in  (2.5). 
If  so,  we  have  A(p,<,x)  =  0  in  (2.9)  (see  (2.7))  for  every  (t,x)  €  and  hence  A(p,t,  x) 
has  0  =  p  —  r  eigenvalues  with  negative  real  part.  In  other  words,  for  all  p. 

As  a  result,  jumps  may  still  occur,  but  only  at  singularities  (impasse  points  in  particular) 
which  of  course  need  not  always  exist  (as,  e.g.,  for  monotone  networks;  see  [DeW]  and  the 
references  therein).  This  case  is  relevant  when  (1.1)  is  replaced  by  (2.1)  for  fixed  p  =  po 
and  (2.1)  already  contains  “all”  the  possible  perturbations  of  (1.1)  (see  Remark  1.1),  in 
other  words,  when  no  physical  constant  is  neglected  in  Note  that  the  problems 
for  p  close  to  po  are  then  perturbations  of  of  the  kind  alluded  to  in  Remark  2.1;  that 
is,  those  that  do  not  involve  any  physical  parameter  not  already  represented  in  □ 

An  important  aspect  of  the  concept  of  P-consistent  solution  is  its  strong  bias  toward 
evolution  in  ‘positive  time.  Indeed,  the  concept  of  'P^-consistent  point  is  not  invariant  under 
time-reversal  because  the  change  of  variable  t  =  T  —  s,x{t)  =  y{T  —  s)  transforms  (2.1) 
into  the  system 

(2.10)  A(p,T-s)J  =  -G(p,T-s,y), 
and  hence  (1.1)  (with  A{t)  =  A.(0,  t),  G(t,  x)  =  <j(0,<,x))  into 

(2.11)  AiT-s)^  =  -G{T-s,y). 

For  (2.11)  the  manifolds  W  and  are  obtained  from  the  corresponding  manifolds  W  and 
for  (1.1)  via  the  variable  change  s  =  T—t.  Such  a  simple  correspondence  does  not  carry 
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over  to  the  sets  of  7^^-consistent  points.  In  fact,  the  operator  in  (2.9)  is  changed 

into  -  s,x)  (because  of  the  change  of  G  into  —G  in  (2.11))  and  hence,  if  (t,a:) 

is  -consistent  for  (1.1),  then  (T  -  t,x)  will  never  be  consistent  with  the  perturbation 
(2.10)  of  (2.11),  unless  p  =  r  in  (2.5). 

For  problems  in  which  t  represents  physically  the  time  variable,  this  is  unimportant 
since  physical  relevance  goes  along  with  evolution  in  positive  time.  But  even  in  problems 
where  the  variable  t  may  increase  or  decrease,  e.g.,  a  load  parameter,  the  above  remark  is 
not  necessarily  a  sign  of  the  inadequacy  of  the  concept  of  "P-consistent  solution.  Rather, 
it  indicates  that  the  same  pertmbations  (2.1)  cannot  be  used  to  describe  the  physical 
phenomena  of  interest  when  t  both  increases  as  well  as  decreases.  This  claim  can  be 
corroborated  for  elementary  one-dimensional  plasticitiy  models  where,  indeed,  some  sign 
change  in  pmt  of  the  equations  must  occur  upon  switching  from  loading  to  unloading  to 
preserve  physical  correctness. 

Remark  2.5:  (Inconsistent  initial  conditions):  The  concept  of  P-consistency  can  be  used 
in  connection  with  the  problems  of  inconsistent  initial  conditions  discussed  in  I,  Section  4. 
In  fact,  the  discussion  there  can  be  repeated  almost  verbatim  when  the  initial  condition 
Xq  is  not  Pp-consistent.  In  tliis  case,  should  be  replaced  by  a  point  Xq  such  that 
X+-X-  e  ker  A(io)  and  is  P^-consistent.  Clearly,  the  -  consistency  requirement 

places  a  further  limitation  upon  the  admissible  values  xj .  Dually,  since  lack  of 
consistency  is  more  frequent  than  lack  of  consistency  in  the  sense  of  Part  I,  there  are  now 
more  points  Xq  from  which  an  initial  jump  should  be  expected.  □ 

We  complete  this  section  with  the  observation  that  P  -  consistency  has  the  desirable 
property  of  being  invariant  under  changes  of  variables  preserving  the  semilinear  structure 
of  the  system  (2.1).  In  other  words,  if  the  system  (2.1)  is  multiplied  by  and  x  is 

changed  into  M{p,t)y  with  smooth  (or  sufficiently  smooth)  €  C?T(R  ), 

thus  becoming  a  transformed  system  P^i,  then  the  P^-consistent  points  (resp.  solutions) 
of  (1.1)  are  converted  into  the  P^-consistent  points  (resp.  solutions)  of  the  unperturbed 
transformed  system  (and  vice-versa).  The  verification  is  straightforward  using  Remark  2.1 
and  will  not  be  given  here. 
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Remark  2.6:  It  is  noteworthy  that  the  eigenvalue  condition  for  the  operator  is 

preserved  even  when  the  operator  does  not  satisfy  the  condition  (2.6).  This 

can  be  used  in  reverse  to  extend  the  validity  of  the  theory  presented  here  to  the  case  when 
(2.6)  is  replaced  by  the  weaker  requirement  that 

kev  C  ker  ^4(0,  t)M(0,<),  x  7, 

for  some  smooth  enough  mapping  M  :  x  J  GT(R").  □ 


3.  Relationship  with  the  Tikhonov-Levinson  Theorem. 


In  order  to  justify  the  definitions  of  the  P-consistent  points  and  of  the  ’P-consistent 
solutions  of  (1.1)  by  making  the  connection  with  the  local  approximation  property  men¬ 
tioned  in  the  Introduction,  we  first  examine  the  particular  case  when  p  =  n  in  (2.5)  and 
A{fj,,t)  has  the  special  form 

(31)  .4^0=1  I"  J, 

where  =  0  and  €  £(R"“'’).  The  system  (2.1)  now  has  the  form 


(3.2) 


=  ?i(A‘,t,x), 
B{^,t)x2  =  ff2(^,t,x). 


The  manifold  W  is  the  zero  set  of  the  mapping  52(0?  •,  Oi  a  point  (i,  x)  is  in  if  and 
only  if 

(3.3)  52(0,f,x)  =  0,  ~(0,t,x)  -  ^(0,t,x)ffi(0,t,x)  €  ige  ^(0,t,x). 

By  I,  Proposition  2.1,  (1.1)  is  geometrically  nonsingular  of  index  one  if  and  only  if  (3.3)  is 
equivalent  to 

(3.4)  ff2(0,  t,  x)  =  0,  -^(0,  t,  x)  €  GZ(R”-"). 

0x2 

Thus,  for  a  given  /z  6  7\  as  in  Definition  2.1,  we  find  that  (f,x)  €  if  (3.4)  holds  and 
the  operator 


(3.5) 


1^(0., X)- 
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e  £(R"-’') 


has  r.  -  r  eigenvalues  with  negative  rea]  part.  If  so,  the  operator  (3.5)  is  invertible  and,  of 
course,  its  inverse 


(3.6) 

has  n  -  r  eigenvalues 


with  negative  real  part. 


(3.7) 


Given  e  >  0,  and  because  of  B{d,t)  —  0,  we  have  —  ^B{e,fJ.,t)  with  B{e,fi,i) 

In  particular,  .5(0, /i, to)  =  |f  (0,to);^,  and  hence  B{e,tx,t)  is  invertible 
for  sufficiently  sraaU  e  >  0  and  t  near  to.  Thus,  for  those  values  of  e  and  t,  the  system 

(3.2)  with  fJ.  replaced  by  efi  is  equivalent  to 

=gi{e(i,t,x) 

[  6X2  =  B(€,/x,t)~^g2(e^,i,x) 

The  key  point  now  is  the  following:  Let  x(i)  denote  the  unique  solution  of  (3.7)  for 
«  =  0  (i.e.  of  (3.2)  with  ^  =  0)  satisfying  a:(to)  =  ^o,  where  (to,xo)  is  such  that 
(3.4)  holds.  Then  the  main  hypothesis  of  the  Tikhonov-Levmson  theorem  (see  e.g.  the 
references  in  [H]  or  [OM]),  which  ensures  that  x(i)  can  be  uniformly  approximated  by 
solutions  of  (3.7)  with  e  >  0  in  some  sufficiently  small  interval  about  to,  is  that  the 
operator  B(0,//,to)-'|^(0,to,xo)  has  n  -  r  eigenvalues  with  negative  real  part.  Since 
^(0^  fo)  =  If  {0,  to)Ai,  this  is  precisely  our  condition  for  the  operator  (3.6).  This  justifies 
the  “approximation”  statement  made  earlier  in  Remark  2.2,  at  least  when  has  the 

form  (3.1)  and  p  —  n  in  (2.5). 


We  now  justify  the  relevance  of  the  7^-consistent  solutions  without  assuming  the  special 
form  (3.1).  But  we  do  continue  (momentarily)  to  require  that  (2.1)  is  an  ODE.  Our 
argument  here  is  simply  that,  locally  in  (p,  t),  reduction  to  the  form  (3.1)  is  always  possible. 
We  skip  a  number  of  straightforward  though  tedious  technicalities  and  focus  on  the  mam 
points.  The  desired  reduction  can  be  accomplished  via  transformations 


U(m,  t)  e  GL(E"),  V{p, t)  e  GI(E”) 

such  that  V{p,i)Aip,t)U{p,i)  has  the  form  (3.1),  followed  by  introducing  the  change  of 
variables  i(f)  =  U{p,t)yit).  In  this  way,  the  right-hand  side  G{p,t,x)  becomes 
V{p,t)[-A{p,t)U{p,t)y  +  Gip,t,Uip,t)y)], 
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and  the  projection  Q{t)  may  be  replaced  by  the  projection  ^(0,  t)Q(t)V’(0,  t)~^.  The 
old  and  new  manifolds  W  and  W'  are  related  via  {t,x)  =  {t,U(0,t)y),  and,  by  using 
Q(<)A(0,<)  =  0  and 


ker  F(0,  t) A(0,  t)U{0,  t)  =  U{0,  t)-^  ker  ^(0,  t), 


(3.8) 

Finally,  Q(t) 


ker  >4(0, r) 

F(0,  t)Q{t)D^GiO,  f ,  17(0,  t)y)U{0,  t\ 


V(0,t)-l  ker>l<0,t)' 


is  transformed  into 


her  4(0, t) 


(3.9)  F(0,  t)Q{t)V{0,  t)-^  (^^(0-  Oa^ 

'  ^  ^  l(/(0,t)-l  ker 4(0,1) 

This  can  be  simplified  by  noting  that  two  of  the  three  terms  arising  in  the  expansion  of 
d{VAU)/dfi  yield  no  contribution.  Indeed,  we  have 

mwwm*)-'  (^(o,t)^)4(o,«)a(o,()|„,.„. . .  =  0, 

because  of  ,j  =  0,  and 


V(0,  t)Q{t)V{0,  t)-^  no,  t)A{Q,  t)  (  f^(0,  t)fi 


because  of  Q{t)A{Q,t)  =  0.  Thus,  (3.9)  reads 

fdA 


=  0, 


ll/(0,e)-^  kerA(0,t) 


and,  together  with  (3.8),  this  shows  that  the  “new”  operator  A(/^,t,  r)  is  simply 

C/(0,t)“^A(/x,t,17(0,i)y)f/(0,t) 


Since  X  =  {/(0,t)?/,  this  operator  has  the  same  eigenvalues  as  A(/x,  t,x)  in  (2.9). 

This  shows  that  the  transformation  x  =  17(0,  t)j/  also  relates  the  old  and  new  V^- 
consistent  points,  and  hence  the  transformation  x{t)  —  17(0,  t)y(<)  relates  the  old  and 
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new  "P^-consistent  solutions.  Since  the  solutions  of  the  old  and  new  perturbed  problems 
are  in  correspondence  through  x(t)  =  L/(p,f)y(t),  the  relationship  between  P^-consistent 
solutions  and  solutions  of  the  perturbed  problems  is  the  same  before  and  after  reduction 
of  (2.1)  to  the  form  (3.1). 

Finally,  when  (2.1)  is  a  (completely  nonsingular)  DAE  for  /i  €  R+,  a  similar  reduction 
procedure  changes  the  problem  into  one  in  which  takes  the  form 

(/,.  0  0 

0  0 
0  0  0 


where  now  5(p,t)  €  £(R'’-0  and  B{0,t)  =  0.  Once  again,  the  transformation  preserves 
P^-consistent  points.  The  transformed  system  (2.1)  has  the  form 


(3.10) 


0  =  922i9;t,^)  e 


where  we  have  set  x  =  (xj,X2)  with  X2  —  (^1211^22)  and  g  —  {91,92)  with  g2  {921,922)- 
Points  of  W  continue  to  be  characterized  by  (3.3),  and  geometric  nonsingularity  still 
amounts  to  equivalence  between  (3.3)  and  (3.4).  Thus,  {dg2/dx2){0,t,x)  is  invertible 
at  each  point  (t,x)  €  W^.  Such  a  point  is  P^-consistent  if  and  only  if  the  operator 


(3.11) 


d92 

dx2 


{0,t,  x) 


'11(0,*)/*  o' 

0  0; 

has  p  -  r  eigenvalues  with  negative  real  part  for  p  €  P,  which,  in  particular,  implies  that 
|J(0,t)p  is  invertible.  The  eigenvalues  of  (3.11)  are  the  numbers  A  6  C  for  which  there  is 
hi  e  C"-’’  \  {0},  hi  =  {hii ,  hii )  such  that 

, /l!5-(o.*.**)  IS(o. \  ^  ((if (0.0/*) '•2* 
is(o.‘,oj  V  0 

Assuming  that  {d922/dx22){0,i,x)  is  invertible,  we  conclude  that  the  nonzero  eigenvalues 
A  are  characterized  by  the  property  that  1/A  is  an  eigenvalue  of  the  operator 

^  dg22 


(3.12) 


dgn  _  dgn  f  0922 
dx2i  0X22  \Ox22 


0X21 


(0,t,x)  e£(M'‘~")- 
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The  condition  that  the  operator  (3.12)  has  p  -  r  eigenvalues  with  negative  real  part  at 
some  point  (<0)^o)  £  exactly  the  (main)  condition  required  in  the  dissertation  [Y] 

for  the  uniform  approximation,  in  the  vicinity  of  fo?  of  th®  solution  a;(f)  of  the  problem 
(3.10)  with  /i  =  0,  satisfying  x{to)  =  xq,  by  solutions  of  (3.10)  with  p  replaced  by  e/i,  e  >  0. 
Thus,  once  again,  our  characterization  of  'P-consistent  points  and  definition  of  P-consistent 
solutions  relates  very  closely  to  (a  generalization  of)  the  Tikhonov-Levinson  theorem. 

Remark  3.1:  Consistent  with  the  nature  of  the  arguments  presented  in  this  section,  it 
is  quite  natural  to  ask  when  a  P^-consistent  solution  of  (1.1)  is  the  limit  of  (classical) 
solutions  of  (2.1)  with  /x  replaced  by  e/i,e  >  0,  in  the  sense  of  distributions.  We  do  not 
have  a  general  answer  to  this  important  question.  Note,  however,  that  the  addition  of 
this  (global)  requirement  for  the  admissibility  of  solutions  of  (1.1)  can  only  confine  further 
the  set  of  P^ -consistent  solutions.  As  observed  in  Section  2,  this  set  is  in  general  already 
countable,  and  often  finite,  when  an  initial  condition  is  prescribed.  This  is  to  say  that, 
at  this  stage,  there  is  only  a  small  number  of  ccindidates  left  that  may,  or  may  not,  pass 
the  aforementioned  global  approximation  test.  (The  difficulty  in  the  global  approximation 
question  is  as  follows:  unlike  the  continuous  solutions  of  P^,  the  discontinuous  solutions  of 
(1.1)  are  insensitive  to  the  driving  vector  field  between  the  initial  and  terminal  points  of 
the  jump.  Thus,  solutions  of  (1.1)  may  sometimes  be  continued  by  jumping  over  attractors 
which  otherwise  capture  the  solutions  of  P^.)  □ 

4.  Examples  from  Nonlinear  Circuit  Theory. 

Example  1:  Our  first  example  is  the  simple,  well-studied  circuit  (see  e.g.  [CD]  and 
[SDe])  of  Figure  4.1 


X2 - 

Xi  =  X23  -  X2 

Figure  4.1 
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which  is  modelled  by  the  autonomous  DAE 

Oil  =  —X2, 

(4-1)  3 

0  =  Xi  —  ^2  +  2:2 . 

The  solutions  are  constrained  to  the  cubic  path  {x  €  :  Xj  =  xf  —  X2}  of  Fig¬ 

ure  4.2.  It  is  trivial  to  check  that  (4.1/2)  is  geometrically  nonsingular  of  index  one. 
In  line  with  I,  Remark  2.4  (change  of  notation  for  autonomous  problems),  the  points 
=  (~2/(3v^),1/\/5)^,  and  ^2  =  “6  are  accessible  impasse  points.  As  discussed, 
e.g.,  in  [CD],  or  [SDe],  at  these  points  the  solution  is  certainly  expected  to  exhibit  jump 
phenomena. 

•1 

0 

-1 

Figure  4.2 

For  the  circuit  of  Figure  4.1  it  is  well-known  (see  loc.  cit.)  that  a  small  parasitic  line 
inductance  //  >  0  in  series  with  the  nonlinear  resistor  should  not  be  disregarded.  The 
resulting  perturbed  circuit  in  Figure  4.3  is  modelled  by  the  DAE 

Cxi  =  -X2, 

(4.2)  ^X2^Xi~X3, 

0  =  X3  -  X2  -I-  X2. 

Evidently,  for  //  =  0  this  becomes  an  augmented  version  of  (4.1)  in  which  the  trivial 
equation  0  =  Xi  —  X3  has  been  added  and  the  variable  Xi  in  the  last  equation  of  (4.1) 
has  been  replaced  by  X3.  It  should  be  evident  that  (4.1)  and  its  augmented  form  are  fully 
equivalent  except  that  the  solutions  are  now  constrained  to  the  cubic  path  W  =  {x  €  R®  : 
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xx  =  X3,  X3  =  X2  ~  ^2}  K®-  The  points  and  ^2  become  the  points  |i  and  ^2  of  W, 

where  =  (-2/3  A  V  A,  -2/3  A)^  and  ^‘2  =  “li- 


Figure  4.3 


Again,  in  line  with  I,  Remark  2.4  it  is  readily  checked  that,  for  /^  >  0,  we  have  W^(fi)  = 
{x  e  :  X3  —  X2  4-  a;2  =  0)  and  that  (4.2)  is  geometrically  nonsingular.  The  operator  A 
of  (2.9)  has  here  the  matrix  representation 


A(M,a:) 


Thus,  for  fi  ^C+  —  (0, 1)  (see  (2.2)),  the  set  of  -consistent  points  is 


W=’>^:^{xeW:\x2\>^}. 


This  excludes  jumps  from  any  point  on  the  upper  or  lower  branches  of  the  5-shaped  curve 
W  except,  of  course,  or  ^2-  From,  say,  the  jump  has  to  end  at  a  point  x"*"  €  kF 
such  that  x"*"  —  |i  =  (0,7?,C)^  G  kerA(O).  Hence  by  (4.2)  with  p  =  0  we  have  C  =  0  and 
therefore  x+  =  (—2/(3  A),  2:2, —2/(3  A))^-  In  other  words,  the  jump  is  parallel  to  the 
X2-axis.  A  similar  conclusion  holds  for  ^2-  The  corresponding  results  for  (4.1)  2ire  obtained 
by  projections  onto  the  plane  X3  =  0  and  shown  in  Figure  4.2. 


These  results  agree  fully  with  those  of  [SDe]  where  the  (standard)  singularly  perturbed 
ODE-form 


Cxi  =  -X2, 

pX2  =  Xi  —  X2  -f  X2, 
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is  used.  This  form  is  obtained  by  elerainating  the  variable  X3  from  the  second  differential 
equation  in  (4.2)  by  using  the  algebraic  equation.. 

As  discussed  in  [SDe],  the  equivalence  of  the  singular  perturbation  (4.3)  with  the  per¬ 
turbation  obtained  by  adding  a  parasitic  inductance  in  series  with  the  nonlinear  resistor 
rests  upon  the  fact  that  this  resistor  is  driven  by  the  current  X2.  An  analogous  equivalence 
pertains  in  circuits  where  a  small  capacitor  is  added  in  parallel  to  a  voltage  driven  resistor 
(see  also  Example  3  below).  But  clearly,  when  the  characteristics  of  a  nonlinear  port  are 
not  globally  solvable  in  terms  of  some  of  the  driving  variables,  then,  as  noted  in  [SDe],  this 
is  symptomatic  of  neglected  parasitics  inside  the  port  whose  effects  cannot  be  reflected  at 
the  terminals.  In  this  case,  we  can  hardly  expect  a  perturbation  by  means  of  the  parasitics 
to  be  equivalent  with  a  singularly  perturbed  ODE.  Furthermore,  no  general  principle  is 
hidden  behind  the  remark  that  (4.3)  can  be  derived  from  (4.1)  by  formaUy  substituting 
for  0  in  the  algebraic  part  of  (4.1).  To  see  this,  note  that  (4.1)  is  unchanged  if  the 
equation  0  =  n  -  +  X2  is  replaced  by  0  =  x^  -  X2  -  while  the  system 

Cxi  =  -X2, 

(4  3')  0  =  x®  -  X2  -  Xi, 

is  one  for  which  the  P,-consistent  points  of  (4.1)  are  now  those  corresponding  to  the  middle 
branch  of  the  S-curve  in  Figure  4.2.  That  the  system  (4.3’)  is  physicaUy  meaningless  is 
clear  from  the  fact  that /iX2  =  x^-X2-xi  amounts  tox^-X2-Xi  =  0  and  /iX2  -X3-X1; 
that  is,  the  current  X2  flows  in  the  direction  opposite  to  physical  reality  (see  Figure  4.1). 

Example  2:  The  first  example  was,  of  course,  an  idealized  circuit  and,  in  practice,  the 
nonlinear,  cubic  resistor  corresponds  to  a  port  made  up  of  several  resistors  and  transistors. 
As  an  example  of  such  a  more  reahstic  version  Example  1  we  consider  the  two-transistor 
port  of  Figure  4.4a  which,  as  before,  is  combined  with  a  capacitor  into  the  circmt  of  Figure 
4.4b.  The  port  was  discussed  in  [CYY]  and  also  (TWj.  We  use  as  independent  variables 
the  voltage  Xi  across  the  port,  the  current  X2  through  the  port  and  the  voltage  differences 
X3  -  -  U4,  X4  =  U3  -  U4,  X5  =  -  Ui,  X5  =  U6  -  U5,  where  U3,  •  •  ■  ,  V6  are  the  voltage 

levels  at  the  indicated  points  of  Figure  4.4. 


20 


Figure  4.4 

Then  the  circuit  is  modelled  by  an  autonomous  DAE  of  the  form 

Cxi  =  X2, 

(4.4) 

0  =  Lx  +  K^(x), 


where 


(4.5) 


Here  the  Ebers-Moll  transistor  model  EMI  was  applied  and  we  set  ^6(s)  =  lEslexp(qs)-l], 
where  I^s  denotes  the  emitter-base  saturation  current,  and  the  constant  q  depends  on  the 
temperature.  For  the  transistor  constants  we  used  aR  ~  0.5,  ap  =  0.99,  Ips  =  10“^®,  and 
q  =  40.  For  simplicity,  in  all  the  calculations  we  worked  with  C  =  1. 


21 


The  DAE  (4.4/5)  has  a  similar  behavior  as  (4.1).  The  Jacobian  of  the  algebraic  part  has 
the  form  L  +  KD^x)  where  D^x)  =  diag  (0, 0,  <;i'(x3), . . .  ,  ^^'(xg)).  An  evaluation  with 
Mathematica  shows  that  when  the  sixth  colum  is  deleted  the  determinant  of  the  resulting 
5x5  matrix  is  a  polynomial  in  the  positive  quantities  (^'(xs),. ..  ,';^'(x6)  which  has  only 
positive  coefficients.  Thus  this  submatrix  is  nonsingular  for  all  x  e  R®  which  implies  that 
the  algebraic  solution  set  is  a  1-dimensional  submanifold  of  R®.  This  path  was  computed 
by  means  of  the  continuation  code  PITCON  ([BRh])  started  from  x  =  0  and  its  projection 
onto  the  Xi,X2'plane  has  the  5-shaped  form  of  Figure  4.5. 


Figure  4.5 

Similarly  it  can  be  shown  that  the  DAE  is  geometrically  nonsingular  of  index  one:  By 
using  the  first  (linear)  algebraic  equations  xi  can  be  we  eliminated  from  the  remaining 
algebraic  equations.  The  resulting  system  has  the  form  TLx  +  TK^{u)  —  0  where  F  is 
some  Gaussian  transformation  matrix.  Now  consider  the  4x4  matrix  obtained  by  deleting 
the  first  row  as  well  as  the  first  and  last  columns.  As  before,  a  Mathematica  evaluation 
shows  that  the  determinant  of  this  matrix  is  a  polynomial  in  ^'(^a),-  •  •  >  which  has 

only  negative  coefficients.  Hence  the  matrix  is  nonsingular  for  all  x  e  R®  and  the  system 
is  nonsingular  of  index  one.  On  the  solution  path  we  observe  an  accessible  as  well  as  an 
inaccessible  singular  point  with  projections  and  ^2i  respectively,  onto  the  Xi,X2-plane. 

As  before,  with  a  small  parasitic  inductance  //  >  0  in  series  with  the  nonlinear  port,  the 
resulting  circuit,  shown  in  Figure  4.6,  is  modelled  by  the  DAE 

(4.6)  Cxi=X2,  ^X2  =  xi-i7,  Ly  +  K^{y)  =  0,  y  =  (x2, . . .  .xt)^ 
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where  L,  K,  $  axe  as  in  (4.5).  Evidently  for  p  —  0  we  obtain  an  augmented  form  of  (4.4) 
which  is  equivalent  with  it  when  the  identification  Xi  =  X7  is  used.  The  projection  of  any 
solution  onto  the  xi,X2-plane  is  again  the  path  of  Figure  4.5. 


Figure  4.6 

For  fi  >  0  the  solutions  are  constrained  to  the  two-dimensional  manifold  W(fi)  =  {x  € 
:  Ly  +  K^{y)  =  0},  y  =  (x2,...  As  in  the  case  of  the  unperturbed  DAE 

(4.4)  it  can  be  shown  that  (4.6)  is  geometrically  nonsingular  of  index  one.  The  matrix 
representation  of  (2.9)  is  here 


(4.7)  A(/x,x)  =  Mdiag  (;/,  0,0, 0,0,0),  M - 


L  +  KD^x] 


,  =  (1,0,0, 0, 0,-1) 


A  computational  determination  of  the  single  nonzero  eigenvalue  during  the  computational 
trace  of  the  manifold  with  PITCON  showed  that,  as  in  the  first  example,  the  set  of  Vfi 
-  consistent  points,  e  =  (0, 1),  consists  of  all  points  projecting  onto  the  upper  and 
lower  branches  of  the  S'-curve  shown  in  Figure  4.5  (excluding,  of  course  the  impasse  points). 
Thus  we  expect  only  a  jump  from  onto  the  (unique)  point  x  on  the  lower  branch  with 
the  same  xi -coordinate.  Figure  4.5  already  shows  an  approximating  solution  for  fi  =  0.01 
computed  for  (4.6)  with  the  IRK-solver  Radau5  of  [HW]  starting  from  a  point  on  the  upper 


A  different,  no  longer  physically  motivated,  perturbation  of  the  original  system  (4.4/5) 
is  given  by  the  corresponding  singularly  perturbed  ODE 


diag  (C,p,...  ,/i)  = 


,  a^  =  (0,1, 0,0, 0,0). 


\Lx  +  K^x)J 

Thus,  in  this  case,  A(^,  x)  is  the  inverse  of  the  matrix  obtained  by  deleting  the  first 
column  of  Lx  -|-  J\D$(x).  The  eigenvalues  of  this  matrix  were  computed  during  a  trace 
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of  the  solution  manifold  and,  now,  it  was  found  that  none  of  the  points  on  the  solution 
manifold  is  -p-consistent.  In  other  words,  this  singular  perturbation  provides  a  different 
(and  erroneous)  information  about  the  admissibility  or  inadmissibility  of  the  solutions  of 
(4.4). 
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Figure  4.7 

Example  3:  Earlier  we  stressed  the  necessity  of  allowing  the  perturbed  system  to 
involve  a  vector-valued  /i.  As  an  illustration  we  consider  the  circuit  of  Figure  4.8  which 
can  be  modelled  by  the  DAE 


(4.8) 


/  -X3  X 

Ai  =  G{x),  A  =  diag(C,  0, 0),  G{x)  =  xi+Xi  +  xs-xl 

\  X3  —  X2  +  X2  / 


Thus,  in  the  terminology  of  [CYY]  we  have  here  an  5-type  and  an  iV-type  nonlinear 
resistors  in  series.  The  solutions  are  constrained  to  a  one- dimensional  submanifold  of 
R*  which  has  the  global  parametrization  x  —  ^^(cr)  =  (7r(cr),<7,  o'  —  <r^)^  where  7r(<T)  = 


-2cr -h  +  (o- -  <7*)^ 
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Figure  4.8  shows  parts  of  this  curve.  The  system  has  three  equilibria,  namely  for  cr  =  0,  ±1. 
Moreover,  there  are  singularities  at  the  points  corresponding  to  zeros  of  the  derivative  of 
TT.  These  Eire  the  four  points  where  =  ±0.854398  . . . ,  cr^  =  ±1.18505 ....  In  Figure  4.8 
the  equilibria  and  singular  points  on  that  part  of  the  curve  are  marked  by  j]i  ,  r?2  and  , 
^2,  respectively,  and  the  arrows  show  the  direction  of  the  flow  between  these  points. 


Now,  in  addition  to  a  parasitic  inductance  in  series  with  the  resistors,  we  introduce, 
as  in  [GW],  also  a  small  capacitance  in  parallel  with  the  iV-type  resistor.  The  resulting 
circuit  of  Figure  4.9  is  modelled  by  the  DAE 


(4.9)  A(fi)i  =  G(x),  A(/i)  =  diag(C,/i,,P2,0,0,0), 


/  \ 


-Xs 


G(x)  = 


Xi  ~Xi 
Xi  ±X2  ±X3 

16  -  X2  ±  xl 

\  X3~  X5-X6  / 


Once  again,  for  ^  this  becomes  an  augmented  version  of  (4.8)  which  is 

equivalent  with  that  DAE  under  the  identification  X4  =  xi,  X5  =  0,  xe  =  xs-  Clearly,  there 
is  little  reason  to  assume  that  both  pi  and  H2  should  be  equal.  Evidently,  the  constraint 
manifold  W(p)  is  here  a  three  -  dimensional  submanifold  of  E®  and  it  is  readily  verified 
that  (4.9)  is  completely  nonsingular. 

A  Mathematica  calculation  shows  that  for  points  ^(cr)  €  W(0),  cr  €  E,  and  with  v  = 
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^2/mi5  the  nonzero  eigenvalues  of  the  operator  (2.9)  are  given  by 

(4.10)  A±  =  -;T-^[g(<T,:/)±((?(a,i/)V4i/7r'((T))i/^],  q(a,u)  =  \  +  v{\-^(x'^)-{(T-a^f. 

ZTT  (a) 

Clearly,  the  "P^-consistent  points  depend  here  on  the  choice  of  //  €  C+  =  (0, 1)^.  Table  4.1 
shows  the  situation  for  i/  =  1;  that  is,  for  fn  =  /X2,  and  cr  <  0  (coresponding  to  the  part  of 
^(0)  in  the  positive  half-plane). 


a-interval 

eigenvalues 

1 

(-0.746664,0) 

complex 

positive  real  parts 

2 

(-0.849913,-0.746664) 

complex 

negative  real  parts 

3 

(-0.854398,-0.849913) 

real 

both  negative 

4 

(-1.18505,-0.854398) 

real 

one  positive,  one  negative 

5 

(-1.33885,-1.18505) 

real 

both  negative 

6 

(-1.46455,-1.33885) 

complex 

negative  real  parts 

7 

(-00,  -1.46455) 

real 

both  negative 

Table  4.1,  i/  =  1 

The  table  shows  that  the  set  ^(0)'=’''  of  -  consistent  points  where  fii  =  M2  consists  of 
the  points  <j>{a)  €  W(0)  with  the  parameter  <7  in  the  intervals  numbered  2,  3,  5,  6,  7.  Note 
that  at  the  common  endpoint  of  intervals  1  and  2  a  pair  of  conjugate  complex  eigenvalues 
passes  through  the  imaginary  axis.  Such  points  were,  by  explicit  assumption,  ruled  out  in 
[SDe].  Moreover,  the  P^  -  consistent  intervals  2,  3  correspond  to  points  of  VF(0)  between 
the  equilibrium  t]i  and  the  singular  point  on  the  curve.  This  shows  that  it  would  be 
fallacious  to  conclude,  say,  on  the  basis  of  the  first  two  examples,  that  for  problems  with  a 
one- dimensional  constraint  manifold  the  boundary  of  the  set  of  P^-consistent  points  always 
coincides  with  the  singular  points  of  the  unperturbed  DAE. 


<T-interval 

eigenvcilues 

1 

(-0.311604,0) 

real 

both  positive 

2 

(-0.593199,-0.311604) 

complex 

positive  real  part 

3 

(-0.722998,  -0.593199) 

complex 

negative  real  part 

4 

(-0.854398,-0.722998) 

real  i 

both  negative 

5 

(-1.18505,-0.854398) 

real 

one  positive,  one  negative 

6 

(-00,-1.18505) 

real 

both  negative 

Table  4.2,  u  =  10.0 
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In  order  to  obtain  a  picture  of  the  change  of  -  consistency  when  v  =  ^2! P-i  changes 
we  give  in  Tables  4.2  and  4.3  the  cases  v  -  10.0  and  v  =  0.001,  respectively  Note  that  the 
■p^-consistent  intervals  3,  4  corresponding  to  the  points  between  rji  and  has  grown.  It  it 
is  easily  checked  that  for  1/  0  the  right  endpoint  of  the  interval  of  "P^-consistent  points 

tends  to  s  =  1  /  \/3- 


<T-interval 

eigenvalues 

1 

(-0.854398,0) 

real 

both  positive 

2 

(-1.18505,-0.854398) 

real 

one  positive,  one  negative 

3 

(-1.20948,-1.18505) 

real 

both  positive 

4 

(-1.21434,-1.20948) 

complex 

positive  real  part 

5 

(-1.22013,  -1.21434) 

complex 

negative  real  part 

6 

(-00,-1.22013) 

real 

both  negative 

Table  4.3,  1/  =  0.001 

Now  the  interval  of  P^-consistent  points  between  tji  and  has  disappeared  and  there 
are  two  small  P^  -  consistent  intervals  3,  4  corresponding  to  points  on  the  upper  branch 
of  the  curve  beyond  the  singular  point  ^2  •  Moreover,  the  common  point  of  intervals  4  and 
5  is  again  a  point  where  a  pair  of  complex  eigenvalues  passes  through  the  imaginary  axis. 


Figure  4.10 


It  is  interesting  to  analyze  the  behavior  of  the  solutions  of  the  perturbed  problem  (4.9) 
near  this  point.  We  expect,  of  course,  a  discontinuity.  But  since  in  this  case  there  is  no 


further  Vft  -  consistent  point  on  VF(0)  in  the  intersection  with  ker  A(0),  there  cannot  be 
any  jump  to  another  branch.  Figure  4.10  shows  results  of  a  computation  with  RadauS 
which  appear  to  suggest  that  the  solution  of  the  perturbed  system  with  /Z]  =  e  =  10“^ 
and  =  i/£,  i'  =  0.001  becomes  periodic.  But  this  is  misleading.  In  fact  the  solution  is 
a  spiral  which  slowly  tends  to  infinity.  This  can  be  seen  well  in  Figure  4.11  for  the  case 
1/  =  10.0  (and  the  same  e).  In  that  case,  interval  6  in  Table  4.2  is  "P^-consistent  and  at  the 
singular  point  ^2  the  solution  becomes  the  spiral  of  the  figure. 


Figure  4.11 

The  shape  of  the  spiral  is  easily  explained  by  setting  =  e,  //2  =  in  (4-9)-  Then  by 
introducing  the  scaled  time  r  =  t/e  the  differential  equations  become 

Cxi  =  -ex3, 

X2  =  —15, 

I/Z3  =  Zi  —  X4. 


For  fixed  z/  and  e  ^  0  this  shows  that  xi  becomes  constant  in  agreement  with  our  theory 
which  requires  that  the  jumps  in  the  discontinuous  solution  are  parallel  to  ker  ^(0).  With 
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=  -y  =  ir(— 1.21463)  we  have  then  the  limiting  DAE 


±2  =  —xs, 
vis  = 

0  =  X4  -\-  X2  X3  —  Xg, 

0  =  X6  —  X2  +  2;2t 
0  =  X3  ~  2:5  —  Xfi, 

with  the  perturbation  parameter  v.  For  v  =  Q  the  solution  manifold  is  one-dimensional 
and  has  the  form  of  the  cubic  xg  =  ~T  ~  2:3  -|-  Xj  when  projected  into  the  Xg,  xg-plane.  It 
is  readily  checked  that  the  upper  and  lower  branches  of  this  cubic  are  "Pir-consistent  under 
the  given  perturbation  with  >  0  and  hence  that  the  solution  should  exhibit  jumps  of  the 
same  type  as  shown  in  Figure  4.2.  This  is  indeed  confirmed  in  Figure  4.12  by  projecting 
the  results  of  Figure  4.10  onto  the  X2,X3-plane  and  superimposing,  at  the  same  time,  the 
cubic  path.  But  note  that  for  small  positive  v  we  no  longer  have  X4  =  7  and  hence  that 
the  cubic  is  slowly  transposed  toward  larger  X4  =  Xj  values;  that  is,  that  the  solution  is 
indeed  a  spiral  tending  to  infinity.  This  shows  that  in  problems  with  higher  dimensional 
perturbation  vectors  the  discontinuity  behavior  may  be  very  different  from  the  straight 

jumps  encountered  in  the  earlier  examples. 
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Figure  4.12 

Example  4:  So  far  all  the  examples  involved  autonomous  DAE’s.  This  is,  of  course, 
readily  changed  if,  in  any  of  these  examples,  the  dynamics  is  not  only  supplied  by  some 
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capacitor  or  inductance,  but  also  by  a  time-dependent  voltage  or  current  source.  A  simple 
example  of  this  type  is  shown  in  Figure  4.13.  It  is  modelled  by  a  DAE  of  the  form  (1.1) 
with 

(Xz~X2  \ 

2:3  -  . 

Xi  X2  X2  / 

where  the  input  voltage  v  is  assumed  to  be  a  smooth  function  of  time.  Clearly,  we  have 
here  a  two-dimensional  constraint  manifold  IF  in  R  x  R®.  Although,  in  contrast  to  the 
current-driven  resistor  of  the  circuit  in  Figure  4.1,  the  resistor  is  now  voltage  driven,  it  is 
evident  that  the  projection  of  W  onto  the  xi,2:2'pl3^6>  has  the  form  of  the  cubic  path  of 
Figure  4.2.  But  the  impasse  points  form  the  submanifolds 

(4.14)  {(t,  x)  G  IF  :  xi  =  :f2/(3V3), r2  =  ±l/Vs]. 


Figure  4.13 

In  analogy  to  Figure  4.6  we  consider  the  perturbed  circuit  of  Figure  4.14  with  a  parasitic 
capacitor  in  parallel  to  the  resistor. 


Figure  4.14 
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This  leads  to  a  DAE  of  the  form  (4.9)  with 


(4.15) 


I  i 


A(fi)  =  diag  (L,  pi,  0, 0, 0) 


G(x)  = 


X3  -  v{t) 
X5  +  X2-  x:l 


\ 


X1-X4  -  X5 


/ 


For  /X  =  0  and  with  the  identification  X4  =  0,  X5  =  xi  this  DAE  is  equivalent  with 
(1.1)/(4.13).  Moreover,  for  fx  >  0,  (4.9/14)  is  easily  checked  to  be  geometrically  nonsingular 
of  index  one.  As  in  the  previous  examples,  the  matrix  representation  of  the  operator  A 
is  readily  calculated  and  shows  that  for  fi  €  C+  —  (0, 1),  the  set  of  7^^  -  consistent  points 
consists  of  those  points  x  G  kF  where  |x2| 


Figure  4.15 

A  solution  of  this  perturbed  DAE  (4.9/14)  was  computed  with  RadauS  for  L  =  1, 
v(t)  =  sin(0.2t),  and  the  starting  point  x  =  (0,1, 0,0,0)^  at  time  f  =  0.  The  result  is 
shown  in  Figure  4.15.  As  expected,  when  the  solution  reaches  a  point  on  the  submanifolds 
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(4.14)  of  impasse  points  it  jumps  between  the  outside  leaves  of  W .  When  projected  onto 
the  xi,X2-plane  the  picture  looks  similar  to  Figure  4.2.  But,  as  Figure  4.15  shows,  in 
1 1^X2, 2:3 -space  the  solution  spirals  several  times  until  it  turns  around  at  a  point  where 
X2  =  X3  and  then  begins  a  new  round  of  spirals,  no  longer  shown  in  the  figure.  This 
indicates  clearly,  that  for  higher  dimensional  algebraic  solution  manifolds  we  may  expect 
a  considerably  more  complex  solution  behavior  than  in  the  one-dimensional  case.  In  fact, 
it  is  hardly  surprising  that  in  problems  of  an  appropriate  form  a  chaotic  regime  may  be 
expected  to  develop  as  observed  by  several  authors  (see  e.g.  [EC]). 
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